Electronic structure interpolation via atomic orbitals 
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We present an efficient scheme for accurate electronic structure interpolations based on the sys- 
tematically improvable optimized atomic orbitals. The atomic orbitals are generated by minimizing 
the spillage value between the atomic basis calculations and the converged plane wave basis cal- 
culations on some coarse fc-point grid. They are then used to calculate the band structure of the 
full Brillouin zone using the linear combination of atomic orbitals (LCAO) algorithms. We find 
that usually 16 - 25 orbitals per atom can give an accuracy of about 10 meV compared to the 
full ab initio calculations. The current scheme has several advantages over the existing interpola- 
tion schemes. The scheme is easy to implement and robust which works equally well for metallic 
systems and systems with complex band structures. Furthermore, the atomic orbitals have much 
better transferability than the Shirley's basis and Wannier functions, which is very useful for the 
perturbation calculations. 

PACS numbers: 71.15.Ap, 71.15.Dx 
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I. INTRODUCTION 

Very often one needs to calculate physical properties 
that need very highly dense fc points to get accurate re- 
sults. For example, Yao et al^ used 2 x 10 6 fc points in 
order to get a converged value of intrinsic anomalous Hall 
conductivity. Direct calculations of the properties via the 
first-principles method are too expensive in these cases. 
An efficient electronic structure interpolation method is 
therefore a powerful tool to reduce the computational 
cost while retaining the accuracy at the ab initio level. 

Maximally Localized Wannier Functions (MLWFs) 2 
have been demonstrated as a powerful tool for interpolat- 
ing electronic structures accurately!^ It has been used to 
calculate the anomalous Hall conductivity, 4 low-field Hall 
conductivity, and magnetic circular dichroism spectrum, 3 
etc. In this approach, the band connectivity, band cross- 
ings and avoided crossings are treated correctly. The 
evaluation of the velocity matrix can also be performed 
analytically. Because of their local character, the ML- 
WFs serve as a kind of exact tight-binding basis of min- 
imal sizes. However, in practice, it is not always easy to 
generate highly localized Wannier functions for metallic 
systems, and systems with complex band structures. 

The electronic structures can also be interpolated by 
the Bloch states. Shirley proposed an interpolation 
scheme in 1996£ The optimal basis set is extracted from 
the full ab initio calculations on a coarse fc-point grid. 
It has been shown that by using about 35 basis func- 
tions per atom, the optimal basis can reproduce the band 
structures of Si and LiF within an accuracy of 10 meV. 
Recently, Shirley's scheme has been extended to more 
general fc-point sampling^ The Shirley's method can be 
used well for both insulators and metals. However, the 
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Shirley's basis is delocalized, therefore it can not take the 
advantages of the locality as MLWFs. Furthermore, the 
Shirley's scheme requires storing a large number of wave 
functions, which might be problematic for large systems. 

Atomic orbitals have been popular recently for self- 
consistent ab initio calculations. Relatively compact ba- 
sis size can give quite satisfactory results compared to the 
plane wave calculations. 7-10 However, to get extremely 
high accuracy band structures (errors ~ 2 - 3 meV), the 
small, rigid atomic basis sets are obviously not adequate. 
Recently the authors (CGH) developed a method to gen- 
erate systematically improvable fully optimized atomic 
orbitals for ab initio calculations »ii The CGH scheme can 
generate orbitals with arbitrary angular momentum and 
multi-zeta functions in a unified procedure. The CGH or- 
bitals offer highly accurate and transferable atomic basis 
sets, which have been tested for a wide variety of sys- 
tems, including semiconductors, oxides, metals, clusters, 
etc. 

We show in this paper that the CGH method can be 
revised to give an efficient scheme for electronic struc- 
tures interpolations. This scheme keeps the advantages of 
both MLWFs and Shirley's methods and gets rid of some 
of their disadvantages. The atomic orbitals are highly 
localized (more localized than MLWFs), and therefore 
can take the advantage of the locality. The atomic or- 
bitals are very neat, and only the one dimensional radial 
functions are stored. We show that the scheme works 
well for both metals and insulators. The highly localized 
orbitals can always be generated for metallic systems. 
Additionally, the transferability of local orbitals for inter- 
polation purposes is explored. Once the atomic orbitals 
are generated, the electronic properties can be calculated 
efficiently by the well established linear combination of 
atomic orbitals (LCAO) algorithms, with the accuracy of 
ab initio plane wave calculations. 

The rest of the paper is organized as follows. In Sec.lHl 
we introduce briefly the method we use to do the band 
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structure interpolation via atomic orbitals. In Sec. lIIIi we 
show, via a few examples, how the present interpolation 
scheme works for both metallic and insulating systems, 
and systems with complex band structures. We also dis- 
cuss the transferability of our interpolation bases. We 
conclude in Sec. IIVI 



II. METHODS 

All the calculations performed in this paper are based 
on the density functional theory (DFT) within local 
(spin) density approximation [L(S)DA]. Although our 
method is demonstrated for the norm-conserving^ pseu- 
dopotentials, in principle, it should also work well for 
the ultrasoft or PAW pseudopotentials. Monkhorst-Pack 
k points^ are used in the following calculation. 



A. Generating local orbitals 

The details of how to construct systematically opti- 
mized atomic basis for ab initio calculations have been 
given in Ref. ITU Here we give a brief review of our ear- 
lier approach to generating local orbitals. The atomic 
orbitals are written as radial functions times spherical 
harmonic functions, i.e., ^(r) = f fl ,i(T)Yi m (f). The ra- 
dial function is taken as the linear combination of Spher- 
ical Bessel functions within certain range r c , i.e., 
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where fj, — {a, i, I, m}, and a,i,£, Z,m are the ele- 
ment type, the index of atom of each element type, the 
multiplicity of the radial functions, the angular momen- 
tum and the magnetic quantum number, respectively. 
ji(qr) is the spherical Bessel function, where q satisfy 
ji(qr c )—0. The coefficients c m are chosen to minimize 
the "spillage" ,±2^ 
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between the Hilbert space spanned by the atomic orbitals 
and the plane wave calculations, for selected reference 
states. In Eq. [21 4 r n (k) is the eigenstate of plane wave 
calculations, and N n and Nk are the number of states of 
interest and the number of k-points in the Brillouin zone. 
P(k) is a projector spanned by the atomic orbitals, i.e., 



P(k)=£|^(k))^(k)(^(k)|, 



(3) 



with 



(r|^(k)) = J2 M* - r„ - R)e lk - (r - +R) , (4) 
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where M (r) is atomic orbitals. r p is the atom cen- 
ter of the fi-th orbital, and R is the lattice vector. 
5 M1 y(k) = {(f> n (k) | <^>„(k)) is the overlap matrix element. 
The minimization of the spillage is achieved by using a 
simulated annealing method. Because the orbitals gener- 
ated from this step usually have unphysical oscillations, 
an additional step is added to smooth out the orbitals by 
minimizing the kinetic energy of each orbital^ 

We have performed extensive tests of this scheme for 
a wide variety of systems, including semiconductors, ox- 
ides, metals, and clusters, etc. The results show that 
the obtained atomic bases are very satisfactory in both 
accuracy and transferability, 11 However, to make high 
accuracy electronic structure interpolation, we need to 
modify the scheme for the task. 

In Ref. ITlL where the aim is to do self-consistent ab ini- 
tio calculations in complex chemical environments, to get 
maximal transferability of the basis sets, a set of dimers 
with different bond lengths are selected as the reference 
systems. The valence states are chosen as the reference 
states. However, in order to interpolate the electronic 
structures with extremely high accuracy for a particular 
system, the system itself is used as the reference system. 
In many cases, the interested states lie only in a given 
energy window (i.e., valence bands as well as conduction 
bands), therefore, we can choose these states (instead of 
all the valence bands) as the references states. 

In the self-consistent LCAO calculations, which are 
usually applied to deal with extremely large systems, it 
is important to keep the basis size compact. Often used 
basis are DZP (double zeta orbitals plus a polar orbital) 
etc. The small basis set is not enough for high accu- 
racy electronic structure interpolations. One must use 
much higher angular momentum and more zeta orbitals. 
In these cases, the CGH scheme has obvious advantages, 
because it can generate the fully optimized orbitals with 
arbitrary angular momentum and any number of radial 
functions for a given angular momentum. 

It has been shown in Ref. [ll| that increasing the or- 
bital radius cutoff r c can reduce the spillage value, and 
therefore improve the quality of the orbitals. However, 
a larger r c also means more computational costs. In the 
self-consistent calculations, one has to balance the com- 
putational costs and the accuracy, and choose a modest 
r c . Yet, for band structure interpolations, the systems 
under study are usually considerably smaller, a larger r c 
can be used. 

Once we obtain the atomic orbitals, the electronic 
structures can be calculated efficiently by using well es- 
tablished LCAO methods. 



III. RESULTS AND DISCUSSIONS 

In this Section, we show via a few examples how the 
present band structure interpolation scheme works. We 
start the tests from simple solid structures like Na, Si and 
Al. Na and Al are metallic while Si is a semiconductor. 
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FIG. 1: (Color online) The RMSE of interpolated bands of 
(a) bcc Na; (b) Si diamond structure; (c) fee Al. The black, 
red, and green columns represent the RMSE of Level 1, 2, 3 
atomic basis sets respectively. 



More complex materials, such as iron based supercon- 
ductor BaFe2As2 are then tested. We further show that 
the atomic bases have much better transferability than 
the the Bloch bases and Wannier bases for the purpose 
of band structure interpolation. 



A. Simple Solids 

We take bcc Sodium as the first example. We use en- 
ergy cut-off E C =5Q Ry. A 8 x 8 x 8 Appoint mesh is used 
as the reference states to generate the local orbitals. We 
choose the lowest 6 bands for interpolation. In order 
to quantify the quality of the interpolation, we calculate 
the root-mean-square errors (RMSE) between the inter- 
polated bands and the plane wave calculations, defined 



RMSE = 
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where N is the number of the k points that the com- 
parison has been taken between e PW the eigenvalue cal- 
culated from plane wave basis and e LCAO the eigenvalue 
calculated from atomic orbitals. 



</)-1x10" 5 



1x10" 



t d erwgy d spillage 

i f energy »--■ f spillage 
► g energy g spillage 




1x10 



6 7 8 9 10 11 12 13 14 15 

r(a.u.) 



.1x10"' 



FIG. 2: (Color online) The spillage values and the averaged 
RMSE as functions of radial cutoff r c for the basis sets, whose 
maximal angular momentum are d, f, g, respectively. 



The RMSE of the 6 interpolated bands of Na are shown 
in Fig.[lja). Three levels of atomic orbitals are generated 
and compared. At each level, we use one group of s, p, 
d, and f orbitals (lslpldlf). The radius cutoff r c is cho- 
sen as 12 a.u. We find that the quality of the basis set is 
much better if all the orbitals of the same level are gener- 
ated simultaneously than those generated separately. As 
shown in Fig. [Ha), the RMSE of the lslpldlf, 2s2p2f2f 
and 3s3p3d3f local basis sets are less than 8 meV, 3 meV 
and 2 meV respectively. 

Figure [5] illustrates how the angular momentums and 
radial cutoff affect the interpolation accuracy. The black 
lines indicate the spillage values while the red lines in- 
dicate the averaged RMSE of the 6 bands. The radius 
cutoff changes from 6 a.u. to 15 a.u.. The largest an- 
gular momentum of the orbitals is set from / = 2 (d 
orbitals) to I — 4 (g orbitals). We use 5 radial orbitals 
(zeta functions) for each angular momentum, which en- 
sures that the accuracy improved by multi-zeta orbitals 
is converged. We find consistent change between spillage 
values and the RMSE, confirming that the spillage is a 
good criterion for the interpolation quality. The RMSE 
reduce considerably as r c increases from 6 a.u. to 9 a.u., 
while only small improvement is found when r c is larger 
than 9 a.u.. The d orbitals are not sufficient for the 
high accuracy interpolation, because the corresponding 
RMSE are around 50 to 100 meV even a very large r c 
is used. Adding in f orbitals improve greatly upon d or- 
bitals. For example, for r c =10 a.u., adding f orbitals to 
the interpolation basis reduces the RMSE to 0.86 meV. 
The g orbitals improve further upon f orbitals, but not 
as much as f orbitals to d orbitals. Since the spillage can 
be as small as 10 -6 , one can get highly accurate physical 
properties related to the wave functions. 

The interpolation results for Si diamond structure and 
fee Al are shown in Fig. QJb)(c) respectively. For these 
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show our electronic structure interpolation scheme also 
works well for materials like BaFe2As2 which have com- 
plex band structures. We interpolate the 14th to 24th 
bands, which span about a 6 eV energy window around 
the Fermi level. Figure [3ja) depicts the RMSE for 
each band using different numbers of numerical orbitals. 
When using the lslpldlflg orbital (25 basis per atom), 
the RMSE of all the bands are under 6 meV except for 
the last one, which is 15.7 meV. Once the second radial 
functions are used, the RMSE of all bands are reduced 
to below 1.3 meV except the last one which is 2.48 meV. 
When using the third radial functions, all the RMSE arc 
under 0.6 meV except the last one which is 1.23 meV. 

Figure [3jb) depicts the band structures of BaFe 2 As2 
along the T-Z-X-r-Z line. The band structures cal- 
culated by plane waves are shown in blue dotted line, 
whereas the band structure calculated by 3s3p3d3f3g nu- 
merical orbitals are plotted in black solid lines. The in- 
terpolated bands and the full ab initio calculations are 
essentially indistinguishable. 

We also test GaAs, fee Cu and CuZn alloys, and the re- 
sults are all similar to the the examples above. For GaAs 
and fee Cu, the RMSE are less than 10 meV by using the 
lslpldlflg orbitals (25 atomic orbitals per atom) only. 
We find that the 3s3p3d3f orbitals (48 orbitals per atom) 
can interpolate the 60 bands around fermi level of CuZn 
alloys with RMSE ~ 2 meV, which are quite similar to 
the results given in Ref. H, in which 42 bases per atom 
were used to get the same accuracy^ 
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FIG. 3: (Color online) (a) The RMSE of interpolated bands 
for BaFe2As2. The black, red, and green columns represent 
the RMSE of Level 1, 2, 3 atomic basis sets respectively, (b) 
Compare the interpolated band structures (blue lines) using 
3s3p3d3f3g basis to the plane wave calculations (red dots). 



two materials, we find g orbitals must be included to get 
highly accurate interpolation results. For Si, the lowest 
8 bands are interpolated. The lslpldlflg orbitals can 
interpolate the band structures with RMSE are about 10 
meV for all the bands. If the basis set is doubled, the 
RMSE for all 8 bands are smaller than 2 meV, showing 
extremely good accuracy. When it further increases to 
the level 3 orbitals, the RMSE can be reduced within 1 
meV. For Al, if the first lslpldlflg orbitals are used, 
the RMSE can be reduced to 10 meV. When we dou- 
ble the basis to 2s2p2d2f2g, and RMSE of all the bands 
fall within 3 meV. The RMSE can be further reduced to 
less than 1 meV, when 3s3p3d3f3g orbitals are used, also 
showing remarkable accuracy. 



B. BaFe 2 As 2 

The recent discovery of superconductivity in doped 
iron arsenide has attracted great attention*^ Here we 



C. TRANSFERBILITY 

Sometimes we need to calculate the physical properties 
under certain perturbations, such as strain, defects, elec- 
tric field, etc. In these transferable basis set is 
highly desirable. Shirley's Bloch interpolation does not 
have any transfcrbility. Even with a small change of the 
system (e.g., a deformation of the unit cell), the Bloch 
basis needs to be regenerated. It has been shown that 
Wanner functions do not have satisfactory transferbility 
either / 7 i 18 because they have too many details. Here 
we demonstrate that the CGH orbitals have very good 
transferbility for the interpolation purpose. 

We take the (3,3) carbon nanotube (CNT) as an exam- 
ple, which is a metallic armchair nanotube. The system 
contains 12 atoms. We interpolate the 20th to 27th bands 
around the Fermi level. We use the energy cutoff _E C =100 
Ry. We first generate the interpolation basis for the CNT 
at its equilibrium structure. We use s, p, d and f orbitals 
as the interpolation basis, with the radius cutoff r c =10 
a.u. The RMSE are shown in Fig. 0Ja). The conduc- 
tion bands are more difficult to interpolate than valence 
bands, but still for Level 2, 3, 4 orbitals, the RMSE arc 
less than 8 meV, 4 meV and 2 meV respectively. The in- 
terpolated band structures are shown in Fig.^d) in black 
solid lines compared to the plane wave results shown in 
the red dotted lines. The two band structures are hardly 
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FIG. 4: (Color online) The RMSE of the interpolated bands 
for the (3,3) carbon nanotube of (a) normal, (b) with 10% 
stretching, and (c) with 10% compression, structures, (d), 
(e), (f) show corresponding interpolated band structures (blue 
lines) using 4s4p4d4f basis compared to the plane wave results 
(red dots). 



distinguishable . 

Using the same interpolation basis, we calculate the 
band structures for (3,3) CNTs which are stretched 10% 
or compressed 10% along the tube and their RMSE. The 
RMSE are given in Fig. |UJb)(c), and the corresponding 
band structures are shown in Fig.@Ie)(f) respectively. As 



we see, the RMSE almost do not change with respect to 
the CNT stretching and compressing, even though the 
band structures themselves change dramatically. These 
results demonstrate that the CGH orbitals have excellent 
transferbility for the interpolation purpose, which is very 
important when perturbation calculations are needed. 



IV. CONCLUSION 

We have presented an efficient scheme for accurate 
electronic structure interpolations based on the system- 
atically improvable optimized atomic orbitals. The cur- 
rent scheme has several advantages over the existing in- 
terpolation schemes. We find usually 16 - 25 orbitals 
per atom can give accuracy about 10 meV compared 
to the full ab initio calculations, similar to that of the 
revised Shirley's scheme. However, unlike the Shirley's 
method, in which a large number of three-dimensional 
wave functions have to be stored in the calculations, the 
present scheme only needs to store the one-dimensional 
radial functions. The atomic orbitals are highly localized, 
therefore the scheme has many good features as the inter- 
polation schemes based on maximally localized Wannier 
functions. Even though the number of atomic orbitals 
is greater than that of Wannier functions, generally the 
atomic orbitals are more localized than Wannier func- 
tions, and therefore have much less neighboring atoms. 
The scheme is easy to implement and robust, working 
equally well for metallic systems and systems with com- 
plex band structures. Furthermore, the atomic orbitals 
have much better transferability than the Shirley's ba- 
sis and Wannier functions, which is very useful for the 
perturbation calculations. 
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